// RungeKutta4 -- a 4th Order Runge Kutta Solver
// Robert McDougal, 29 March 2007

import java.util.Vector;

public class RungeKutta4 extends DifferentialEquationSolver {
    // pass constructor arguments to general case
    public RungeKutta4(IOdeRhs yPrime, double[] y0, double t0) {
        super(yPrime, y0, t0);
    }
    
    // Algorithm as presented in Nagle, Saff and Snider
    public double[] step(double h) {
        double[] k1 = multiply(h, yPrime.eval(t0, y0));
        double tHalf = t0 + h * .5;
        double[] k2 = multiply(h, yPrime.eval(tHalf, add(y0, multiply(.5, k1))));
        double[] k3 = multiply(h, yPrime.eval(tHalf, add(y0, multiply(.5, k2))));
        double[] k4 = multiply(h, yPrime.eval(t0 + h, add(y0, k3)));
        double[] dy = multiply(1./6, add(k1, add(multiply(2, add(k2, k3)), k4)));
        t0 += h;
        y0 = add(y0, dy);
        return y0.clone();
    }
    
}
